﻿using System;
using System.Text;

namespace SoftwareConsulting.BI11.NumericalMethods.MatrixAlgebra
{
    /// <summary>
    /// Класс СЛАУ в алгебре матриц
    /// </summary>
    public class MatrixAsLinearEquationsSystem
    {        
        private double[,] _rows;        
        private Vector[] _solutions;

        /// <summary>
        /// Конструктор класса
        /// </summary>
        /// <param name="M">коэффициенты левой части СЛАУ (матрица)</param>
        /// <param name="C">коэффициенты правой части СЛАУ (матрица)</param>
        public MatrixAsLinearEquationsSystem(double[,] M, double[,] C)
        {
            int n = M.GetLength(0);
            if (M.GetLength(1) != n)
                throw new BIMatrixNotSquareException(n,M.GetLength(1), null);
            if (C.GetLength(1) != n)
                throw new BIMatrixCantBuildLinearEquationsException(n, C.GetLength(1), null);
            _rows = new double[n, n + C.GetLength(0)];
            for (int i = 0; i < n; i++)
            {
                for (int j = 0; j < n; j++)
                    _rows[i, j] = M[i, j];
                for (int j = 0; j < C.GetLength(0); j++)
                    _rows[i, n + j] = C[j, i];
            }
        }

        /*не CLS-комплиянт/// <summary>
        /// Конструктор класса
        /// </summary>
        /// <param name="M">коэффициенты левой части СЛАУ (матрица)</param>
        /// <param name="C">коэффициенты правой части СЛАУ (вектор)</param>
        public MatrixAsLinearEquationsSystem(double[,] M, double[] C)
        {
            int n = M.GetLength(0);
            if (M.GetLength(1) != n)
                throw new BIMatrixNotSquareException(n,M.GetLength(1), null);
            if (C.Length != n)
                throw new BIMatrixCantBuildLinearEquationsException(n, C.Length, null);
            _rows = new double[n, n + 1];
            for (int i = 0; i < n; i++)
            {
                for (int j = 0; j < n; j++)
                    _rows[i, j] = M[i, j];
                _rows[i, n] = C[i];
            }
        }*/

        /// <summary>
        /// Конструктор класса
        /// </summary>
        /// <param name="M">матрица левой части СЛАУ</param>
        /// <param name="C">вектор правой части СЛАУ</param>
        public MatrixAsLinearEquationsSystem(Matrix M, Vector C) /*: 
            this( M.Components, C.Components)*/
        {            
            double[,] Mc = M.Components;
            double[] Cc = C.Components;
            int n = Mc.GetLength(0);
            if (Mc.GetLength(1) != n)
                throw new BIMatrixNotSquareException(n, Mc.GetLength(1), null);
            if (Cc.Length != n)
                throw new BIMatrixCantBuildLinearEquationsException(n, Cc.Length, null);
            _rows = new double[n, n + 1];
            for (int i = 0; i < n; i++)
            {
                for (int j = 0; j < n; j++)
                    _rows[i, j] = Mc[i, j];
                _rows[i, n] = Cc[i];
            }
        }

        /// <summary>
        /// Находит решение для константного вектора P, применяя обратную подстановку
        /// </summary>
        /// <param name="P"></param>
        private void SolveUsingBackSubstitution(int P)
        {
            int n = _rows.GetLength(0);
            double[] result = new double[n];            
            for (int i = n - 1; i >= 0; i--)
            {
                double x = _rows[i, n + P];
                for (int j = i + 1; j < n; j++)
                    x -= result[j] * _rows[i, j];
                result[i] = x / _rows[i, i];
            }
            _solutions[P] = new Vector(result);
        }

        /// <summary>
        /// Возвращает позицию наибольшей точки вращения для указанного шага вращения
        /// </summary>
        /// <param name="pivotingStep"></param>
        /// <returns></returns>
        private int GetLargestPivot(int pivotingStep)
        {
            double pivot = Math.Abs(_rows[pivotingStep, pivotingStep]);
            int result = pivotingStep;
            double x;
            for (int i = pivotingStep + 1; i < _rows.GetLength(0); i++)
            {
                x = Math.Abs(_rows[i, pivotingStep]);
                if (x > pivot)
                {
                    result = i;
                    pivot = x;
                }
            }
            return result;
        }

        /// <summary>
        /// Выполняет вращение в указанной точке
        /// </summary>
        /// <param name="pivotLocation"></param>
        private void Pivot(int pivotLocation)
        {
            double inversePivot = 1 / _rows[pivotLocation, pivotLocation];            
            int n = _rows.GetLength(0);
            int m = _rows.GetLength(1);
            for (int i = pivotLocation + 1; i < n; i++)
            {
                double r = inversePivot * _rows[i, pivotLocation];
                for (int j = pivotLocation; j < m; j++)
                    _rows[i, j] -= _rows[pivotLocation, j] * r;
            }
        }

        /// <summary>
        /// Выполняет оптимальное вращение в указанной точке
        /// </summary>
        /// <param name="pivotLocation"></param>
        private void PivotingStep(int pivotLocation)
        {
            SwapRows(pivotLocation, GetLargestPivot(pivotLocation));
            Pivot(pivotLocation);
        }

        /// <summary>
        /// Возвращает решение для первого константного вектора
        /// </summary>
        /// <returns></returns>
        public Vector GetSolution()
        {
            return GetSolution(0);
        }

        /// <summary>
        /// Возвращает решение для указанного константного вектора
        /// </summary>
        /// <param name="constantsIndex"></param>
        /// <returns></returns>
        public Vector GetSolution(int constantsIndex)
        {
            if (_solutions == null)
                SolveUsingPivot();
            if (_solutions[constantsIndex] == null)
                SolveUsingBackSubstitution(constantsIndex);
            return _solutions[constantsIndex];
        }

        /// <summary>
        /// Решает СЛАУ, вращая матрицу
        /// </summary>
        private void SolveUsingPivot()
        {
            int n = _rows.GetLength(1);
            for (int i = 0; i < n; i++)
                PivotingStep(i);
            _solutions = new Vector[_rows.GetLength(1) - n];
        }

        /// <summary>
        /// Меняет местами строки P и Q
        /// </summary>
        /// <param name="p"></param>
        /// <param name="q"></param>
        private void SwapRows(int p, int q)
        {
            if (p == q)
                return;            
            int m = _rows.GetLength(1);
            for (int j = 0; j < m; j++)
            {
                double temp = _rows[p, j];
                _rows[p, j] = _rows[q, j];
                _rows[q, j] = temp;
            }            
        }

        /// <summary>
        /// Возвращает как строку
        /// </summary>
        /// <returns></returns>
        public override string ToString()
        {
            StringBuilder sb = new StringBuilder();
            char[] separator = { '[', ' ' };
            int n = _rows.GetLength(0);
            int m = _rows.GetLength(1);
            for (int i = 0; i < n; i++)
            {
                separator[0] = '(';
                for (int j = 0; j < n; j++)
                {
                    sb.Append(separator);
                    sb.Append(_rows[i, j]);
                    if (j == 0) 
                        separator[0] = ',';
                }
                separator[0] = ':';
                for (int j = n; j < m; j++)
                {
                    sb.Append(separator);
                    sb.Append(_rows[i, j]);
                    if (j == n) 
                        separator[0] = ',';
                }
                sb.Append(')');
                sb.Append('\n');
            }
            return sb.ToString();
        }
    }
}
